Load data

load(file = "Rdata/frap_2020_dis_Rx.Rdata")
load(file = "Rdata/CalMapper_activities_fire_2020.Rdata")
load(file = "Rdata/frap_2020_Rx.Rdata")
load(file = "Rdata/NF.Rdata")

It’s too impossible to actually find overlaps of activity acres so I’m going to stick with dissolving them and finding footprint acreage overlap

Dissolve

act_20_diss <- act_2020 %>% 
  st_union()

Transform to match

act_20_diss <- st_transform(act_20_diss, st_crs(frap_20))

Overlap

overlap <- st_intersection(frap_diss_2020, act_20_diss)
overlap_ha <- as.numeric(st_area(overlap)/10000) %>% 
  round()
overlap_ha
## [1] 6219

Find areas of each

frap_20_ha <- as.numeric(st_area(frap_diss_2020)/10000) %>% 
  round()
calmapper_20_ha <- as.numeric(st_area(act_20_diss)/10000) %>% 
  round()
paste("total area of FRAP dissolved polygons is", frap_20_ha, "ha")
## [1] "total area of FRAP dissolved polygons is 12100 ha"
paste("total area of CAL MAPPER dissolved polygons is", calmapper_20_ha, "ha")
## [1] "total area of CAL MAPPER dissolved polygons is 17478 ha"

Plot with National Forests to make sure it looks right

mapview::mapview(list(act_20_diss, NF, frap_diss_2020),
                 col.regions=list("red","blue", "green"),
                 col=list("red","blue", "green"),
                 alpha.regions = 0.4)

Function to make Venn diagram

venn_calmap <- function(){
  grid.newpage()
  draw.pairwise.venn(area1 = frap_20_ha,
                   area2 = calmapper_20_ha,
                   cross.area = overlap_ha,
                   #category = c("", ""),
                   fill = c("red", "green"),
                   cex = 1)
}
venn_calmap()

## (polygon[GRID.polygon.11], polygon[GRID.polygon.12], polygon[GRID.polygon.13], polygon[GRID.polygon.14], text[GRID.text.15], text[GRID.text.16], text[GRID.text.17], text[GRID.text.18], text[GRID.text.19])
file <- "Venn_FRAP_CalMAPPER.jpeg"
grid.newpage()
png(file = file, width = 1000, height = 1000)
venn_calmap()
## (polygon[GRID.polygon.20], polygon[GRID.polygon.21], polygon[GRID.polygon.22], polygon[GRID.polygon.23], text[GRID.text.24], text[GRID.text.25], text[GRID.text.26], text[GRID.text.27], text[GRID.text.28])
dev.off()
## png 
##   2

Examine treatments in FRAP that don’t overlap

frap_only <- st_difference(frap_20, act_20_diss)
## Warning: attribute variables are assumed to be spatially constant throughout
## all geometries
frap_only %>% head()
## Simple feature collection with 6 features and 18 fields
## Geometry type: GEOMETRY
## Dimension:     XY
## Bounding box:  xmin: -69856.8 ymin: 40323.23 xmax: -45281.23 ymax: 81724.03
## Projected CRS: NAD83 / California Albers
##   YEAR_ STATE AGENCY UNIT_ID TREATMENT_             TREATMENT1 START_DATE
## 1  2020    CA    CDF     AEU      10509  Jan 31 2020 Broadcast 2020-01-31
## 2  2020    CA    CDF     AEU      10630     Feb 2020 Broadcast 2020-02-13
## 3  2020    CA    CDF     AEU      10639     Feb 2020 Broadcast 2020-02-11
## 4  2020    CA    CDF     AEU      10780   Mar 3 2020 Broadcast 2020-03-03
## 5  2020    CA    CDF     AEU      10782 March 5 2020 Broadcast 2020-03-05
## 6  2020    CA    CDF     AEU      11060   2020 Fuels Reduction 2020-06-09
##     END_DATE TREATED_AC GIS_ACRES RX_CONSUM PRE_CON_CL POST_CON_C
## 1 2020-01-31       15.0  15.04700         0          0          0
## 2 2020-02-18       22.4  38.83980         0          0          0
## 3 2020-02-20       75.5  75.49810         0          0          0
## 4 2020-03-03       61.8  61.77150         0          0          0
## 5 2020-03-05       16.5  16.47660         0          0          0
## 6 2020-06-09        4.1   4.14155         0          0          0
##   TREATMENT_TYPE Shape_Leng Shape_Area treatment_type Id
## 1              1  1365.8875   60893.22 Broadcast Burn  0
## 2              1  3918.2200  157179.23 Broadcast Burn  0
## 3              1  3918.8345  305529.98 Broadcast Burn  0
## 4              1  2769.0343  249980.46 Broadcast Burn  0
## 5              1  1452.1062   66678.60 Broadcast Burn  0
## 6              1   814.7527   16760.24 Broadcast Burn  0
##                         geometry
## 1 POLYGON ((-45282.52 81663.6...
## 2 MULTIPOLYGON (((-49811.54 5...
## 3 MULTIPOLYGON (((-47047.56 8...
## 4 MULTIPOLYGON (((-47588.27 8...
## 5 MULTIPOLYGON (((-49808.72 5...
## 6 MULTIPOLYGON (((-69740.64 4...
mapview::mapview(list(frap_20, frap_only), col.regions=list("red","blue"),col=list("red","blue"))
frap_only <- frap_only %>% 
  mutate(area_m2 = st_area(frap_only)) %>% 
  arrange(desc(area_m2))
head(frap_only)
## Simple feature collection with 6 features and 19 fields
## Geometry type: GEOMETRY
## Dimension:     XY
## Bounding box:  xmin: -249552.8 ymin: -460799.2 xmax: 243175.9 ymax: 440038.5
## Projected CRS: NAD83 / California Albers
##   YEAR_ STATE AGENCY UNIT_ID TREATMENT_                       TREATMENT1
## 1  2020    CA    CDF     TCU       9655       2020-Winton Broadcast Burn
## 2  2020    CA    PVT     LNU      10644                      Lake Assist
## 3  2020    CA    CDF     RRU       8359 Lake Matthews FMU Broadcast Burn
## 4  2020    CA    OTH     SLU      10956           2020 C234 CAMP ROBERTS
## 5  2020    CA    FWS     SKU       9579                 Lower Klamath RX
## 6  2020    CA    CDF     TCU      10843               2020-New Hogan VMP
##   START_DATE   END_DATE TREATED_AC GIS_ACRES RX_CONSUM PRE_CON_CL POST_CON_C
## 1 2020-02-01 2020-12-31     123.08  6311.960         0          0          0
## 2 2020-02-21 2020-04-16     200.00  4003.440         0          0          0
## 3 2020-05-21 2020-05-21      81.00   929.511         0          0          0
## 4 2020-05-11 2020-05-15     450.00  1208.190         0          0          0
## 5 2020-11-03 2020-11-03     235.00  2452.790         0          0          0
## 6 2020-04-01 2020-04-29      35.81   537.723         0          0          0
##   TREATMENT_TYPE Shape_Leng Shape_Area treatment_type Id        area_m2
## 1              1   50593.42   25543584 Broadcast Burn  0 23350414 [m^2]
## 2              1   25184.50   16201344 Broadcast Burn  0 14570273 [m^2]
## 3              1   10303.44    3761597 Broadcast Burn  0  3432334 [m^2]
## 4              1   14891.05    4889392 Broadcast Burn  0  3079165 [m^2]
## 5              1   17711.31    9926099 Broadcast Burn  0  2366054 [m^2]
## 6              1   10417.09    2176087 Broadcast Burn  0  1614601 [m^2]
##                         geometry
## 1 MULTIPOLYGON (((-34700.92 4...
## 2 MULTIPOLYGON (((-247635.4 1...
## 3 POLYGON ((242853.4 -460283....
## 4 MULTIPOLYGON (((-68366.27 -...
## 5 MULTIPOLYGON (((-141630.4 4...
## 6 MULTIPOLYGON (((-71642.56 1...
write.csv(frap_only %>% st_drop_geometry(), file = "FRAP_not_CalMapper.csv", row.names = F)

Make a Venn diagram using treated_acres columns instead of GIS acreage

frap_20 and act_2020

frap_only_ID <- frap_20 %>% 
  filter(!TREATMENT_ %in% act_2020$TREATMENT_ID)
frap_only_ID %>% summarize(sum(TREATED_AC))
## Simple feature collection with 1 feature and 1 field
## Geometry type: MULTIPOLYGON
## Dimension:     XY
## Bounding box:  xmin: -330662.3 ymin: -493555.5 xmax: 387954.7 ymax: 430452.3
## Projected CRS: NAD83 / California Albers
##   sum(TREATED_AC)                       geometry
## 1          2902.6 MULTIPOLYGON (((-69002.94 -...
overlap_id_ac <- frap_20 %>% 
  st_drop_geometry() %>% 
  filter(TREATMENT_ %in% act_2020$TREATMENT_ID) %>% 
  summarize(sum(TREATED_AC))
overlap_id_ac
##   sum(TREATED_AC)
## 1        11110.57
frap_id_ac <- frap_20 %>% 
  st_drop_geometry() %>% 
  summarize(sum(TREATED_AC))
frap_id_ac
##   sum(TREATED_AC)
## 1        14013.17
calmapper_id_ac <- act_2020 %>% 
  st_drop_geometry() %>% 
  summarize(sum(TREATED_ACRES))
calmapper_id_ac
##   sum(TREATED_ACRES)
## 1           16390.18
venn_calmap_IDs <- function(){
  grid.newpage()
  draw.pairwise.venn(area2 = frap_id_ac,
                   area1 = calmapper_id_ac,
                   cross.area = overlap_id_ac,
                   #category = c("", ""),
                   fill = c( "green", "red"),
                   cex = 0,
                   ext.text = T,
                   euler.d = T,
                   scaled = T)
}
venn_calmap_IDs()

## (polygon[GRID.polygon.29], polygon[GRID.polygon.30], polygon[GRID.polygon.31], polygon[GRID.polygon.32], text[GRID.text.33], text[GRID.text.34], text[GRID.text.35], text[GRID.text.36], text[GRID.text.37])
file <- "Venn_FRAP_CalMAPPER_IDs.jpeg"
grid.newpage()
png(file = file, width = 1000, height = 1000)
venn_calmap_IDs()
## (polygon[GRID.polygon.38], polygon[GRID.polygon.39], polygon[GRID.polygon.40], polygon[GRID.polygon.41], text[GRID.text.42], text[GRID.text.43], text[GRID.text.44], text[GRID.text.45], text[GRID.text.46])
dev.off()
## png 
##   2